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1. Introduction 

The goal of this note is to discuss some of the applications of discrete complex 
analysis to problems in probability and statistical physics. It is not an exhaustive 
survey, and it lacks many references. Forgoing completeness, we try to give a taste 
of the subject through examples, concentrating on a few of our recent papers with 
Dmitry Chelka k, Hugo Duminil-Copin and Clement Hongler [CS081 ICSTM ICSTOl 
IDCS1 0, HS10 . There are certainly other interesting developments in discrete 
complex analysis, and it would be a worthy goal to write an extensive exposition 
with an all-encompassing bibliography, which we do not attempt here for lack of 
space. 

Complex analysis (we restrict ourselves to the case of one complex or equiv- 
alcntly two real dimensions) studies analytic functions on (sub-domains of) the 
complex plane, or more generally analytic structures on two dimensional mani- 
folds. Several things are special about the (real) dimension two, and we won't 
discuss an interesting and often debated question, why exactly complex analysis is 
so nice and elegant. In particular, several definitions lead to identical class of ana- 
lytic functions, and historically different adjectives (regular, analytic, holomorphic, 
monogenic) were used, depending on the context. For example, an analytic func- 
tion has a local power series expansion around every point, while a holomorphic 
function has a complex derivative at every point. Equivalence of these definitions 
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is a major theorem in complex analysis, and there are many other equivalent def- 
initions in terms of Cauchy-Riemann equations, contour integrals, primitive func- 
tions, hydrodynamic interpretation, etc. Holomorphic functions have many nice 
properties, and hundreds of books were devoted to their study. 

Consider now a discretized version of the complex plane: some graph embed- 
ded into it, say a square or triangular lattice (more generally one can speak of 
discretizations of Riemann surfaces) . Can one define analytic functions on such a 
graph? Some of the definitions do not admit a straightforward discretization: e.g. 
local power series expansions do not make sense on a lattice, so we cannot really 
speak of discrete analyticity. On the other hand, as soon as we define discrete 
derivatives, we can ask for the holomorphicity condition. Thus it is philosophically 
more correct to speak of discrete holomorphic, rather than discrete analytic func- 
tions. We will use the term preholomorphic introduced by Ferrand |Fer44] . as we 
prefer it to the term monodiffric used by Isaacs in the original papers |Isa411 IIsa52) 
(a play on the term monogenic used by Cauchy for continuous analytic functions) . 

Though the preholomorphic functions are easy to define, there is a lack of 
expository literature about them. We see two main reasons: firstly, there is no 
canonical preholomorphicity definition, and one can argue which of the competing 
approaches is better (the answer probably depends on potential applications) . Sec- 
ondly, it is straightforward to transfer to the discrete case beginnings of the usual 
complex analysis (a nice topic for an undergraduate research project), but the easy 
life ends when it becomes necessary to multiply preholomorphic functions. There 
is no easy and natural way to proceed and the difficulty is addressed depending on 
the problem at hand. 

As there seems to be no canonical discretization of the complex analysis, we 
would rather adopt a utilitarian approach, working with definitions corresponding 
to interesting objects of probabilistic origin, and allowing for a passage to the 
scaling limit. We want to emphasize, that we are concerned with the following 
triplets: 

1. A planar graph, 

2. Its embedding into the complex plane, 

3. Discrete Cauchy-Riemann equations. 

We are interested in triplets such that the discrete complex analysis approximates 
the continuous one. Note that one can start with only a few elements of the triplet, 
which gives some freedom. For example, given an embedded graph, one can ask 
which discrete difference equations have solutions close to holomorphic functions. 
Or, given a planar graph and a notion of preholomorphicity, one can look for an 
appropriate embedding. 

The ultimate goal is to find lattice models of statistical physics with preholo- 
morphic observables. Since those observables would approximate holomorphic 
functions, some information about the original model could be subsequently de- 
duced. 
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Below we start with several possible definitions of the preholomorphic functions 
along with historical remarks. Then we discuss some of their recent applications 
in probability and statistical physics. 



2. Discrete holomorphic functions 

For a given planar graph, there are several ways to define preholomorphic functions, 
and it is not always clear which way is preferable. A much better known class is 
that of discrete harmonic (or preharmonic) functions, which can be defined on 
any graph (not necessarily planar), and also in more than one way. However, one 
definition stands out as the simplest: a function on the vertices of graph is said to 
be preharmonic at a vertex v, if its discrete Laplacian vanishes: 

= AH(u) := Yl (H(v)-H(u)) . (1) 

v: neighbor of u 

More generally, one can put weights on the edges, which would amount to taking 
different resistances in the electric interpretation below. Preharmonic functions 
on planar graphs are closely related to discrete holomorphicity: for example, their 
gradients defined on the oriented edges by 

F(uv) := H(v) - H(u) , (2) 

are preholomorphic. Note that the edge function above is antisymmetric, i.e. 
F(uv) = —F(vu). 

Both classes with the definitions as above are implicit already in the 1847 work 
of Kirchhoff Kir47 , who interpreted a function defined on oriented edges as an 
electric current flowing through the graph. If we assume that all edges have unit 
resistance, than the sum of currents flowing from a vertex is zero by the first 
Kirchhoff law: 

E f w = °> ( 3 ) 

u: neighbor of v 

and the sum of the currents around any oriented closed contour 7 (for the planar 
graphs it is sufficient to consider contours around faces) face is zero by the second 
Kirchhoff law: 

The two laws are equivalent to saying that F is given by the gradient of a potential 
function H as in ([2]), and the latter function is preharmonic ((TJ) . One can equiv- 
alently think of a hydrodynamic interpretation, with F representing the flow of 
liquid. Then conditions ^ and @ mean that the flow is divergence- and curl-free 
correspondingly. Note that in the continuous setting similarly defined gradients of 
harmonic functions on planar domains coincide up to complex conjugation with 
holomorphic functions. And in higher dimensions harmonic gradients were pro- 
posed as one of their possible generalizations. 
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There are many other ways to introduce discrete structures on graphs, which 
can be developed in parallel to the usual complex analysis. We have in mind 
mostly such discretizations that restrictions of holomorphic (or harmonic) functions 
become approximately preholomorphic (or preharmonic) . Thus we speak about 
graphs embedded into the complex plane or a Riemann surface, and the choice 
of embedding plays an important role. Moreover, the applications we are after 
require passages to the scaling limit (as mesh of the lattice tends to zero), so we 
want to deal with discrete structures which converge to the usual complex analysis 
as we take finer and finer graphs. 

Preharmonic functions satisfying ([T]) on the square lattices with decreasing 
mesh fit well into this philosophy, and were studied in a number of papers in early 
twentieth century (see e.g. [PW23, Bo u261 ILus26| h culminating in the seminal 
work of Courant, Friedrichs and Lewy. It was shown in [CFL28] that solution to 
the Dirichlct problem for a discretization of an elliptic operator converges to the 
solution of the analogous continuous problem as the mesh of the lattice tends to 
zero. In particular, a preharmonic function with given boundary values converges 
in the scaling limit to a harmonic function with the same boundary values in a 
rather strong sense, including convergence of all partial derivatives. 

Preholomorphic functions distinctively appeared for the first time in the papers 
[Isa41|ITs~a52| of Isaacs, where he proposed two definitions (and called such functions 
"monodiffric" ) . A few papers of his and others followed, studying the first definition 
which is asymmetric on the square lattice. More recently the first definition 
was studied by Dynnikov and Novikov [DN03] in the triangular lattice context, 
where it becomes symmetric (the triangular lattice is obtained from the square 
lattice by adding all the diagonals in one direction). 

The second, symmetric, definition was reintroduced by Ferrand, who also dis- 
cussed the passage to the scaling limit |Fer44[ ILF55] . This was followed by exten- 
sive studies of Duffin and others, starting with [Duf56 . 

Both definitions ask for a discrete version of the Cauchy-Riemann equations 
di a F = id a F or equivalently that ^-derivative is independent of direction. Con- 
sider a subregion J7 e of the mesh e square lattice eZ 2 C C and define a function on 
its vertices. Isaacs proposed the following two definitions, replacing the derivatives 
by discrete differences. His "monodiffric functions of the first kind" are required 
to satisfy inside i7 c the following identity: 

F(z + ie)-F(z) = i(F(z + e) - F(z)) , (5) 

which can be rewritten as 

F(z + ie) - F(z) _ F(z + e) - F(z) 
(z + ie) — z {z + e) — z 

We will be working with his second definition, which is more symmetric and also 
appears naturally in probabilistic context (but otherwise the theories based on 
two definitions are almost the same). We say that a function is preholomorphic, if 
inside J7 e it satisfies the following identity, illustrated in Figure Q] 



F(z + ie) - F(z + e) = i (F{z + e(l + 1)) - F(z)) , 



(G) 
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Figure 1. The first and the second Isaacs' definitions of discrete holomorphic functions: 
multiplied by i difference along the vector a is equal to the difference along the rotated 
vector ia. Note that the second definition (on the right) is symmetric with respect to 
lattice rotations, while the first one is not. 



which can also be rewritten as 

F(z + ie)-F(z + e) _ F{z + e(l + i)) - F(z) 
(z + ie)-(z + e) ~ (z + e(l + i))-z 

It is easy to see that restrictions of continuous holomorphic functions to the mesh 
e square lattice satisfy this identity up to 0(e 3 ). Note also that if we color the 
lattice in the chess-board fashion, the complex identity (JSJ can be written as two 
real identities (its real and imaginary parts), one involving the real part of F at 
black vertices and the imaginary part of F at white vertices, the other one - vice 
versa. So unless we have special boundary conditions, F splits into two "demi- 
functions" (real at white and imaginary at black vs. imaginary at black and real 
at white vertices), and some prefer to consider just one of those, i.e. ask F to be 
purely real at black vertices and purely imaginary at white ones. 

The theory of so defined preholomorphic functions starts much like the usual 
complex analysis. It is easy to check, that for preholomorphic functions sums 
are also preholomorphic, discrete contour integrals vanish, primitive (in a simply- 
connected domain) and derivative are well-defined and are preholomorphic func- 
tions on the dual square lattice, real and imaginary parts are preharmonic on their 
respective black and white sublattices, etc. Unfortunately, the product of two 
preholomorphic functions is no longer preholomorphic: e.g., while restrictions of 
1, z, and z 2 to the square lattice are preholomorphic, the higher powers are only 
approximately so. 

Situation with other possible definitions is similar, with much of the linear 
complex analysis being easy to reproduce, and problems appearing when one has 
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to multiply preholomorphic functions. Pointwisc multiplication cannot be con- 
sistently defined, and though one can introduce convolution-type multiplication, 
the possible constructions are non-local and cumbersome. Sometimes, for different 
graphs and definitions, problems appear even earlier, with the first derivative not 
being preholomorphic. 

Our main reason for choosing the definition ([6]) is that it naturally appears in 
probabilistic context. It was also noticed by Duffin that (jSJ) nicely generalizes to a 
larger family of rhombic lattices, where all the faces are rhombi. Equivalently, one 
can speak of isoradial graphs, where all faces are inscribed into circles of the same 
radius — an isoradial graph together with its dual forms a rhombic lattice. 

There are two main reasons to study this particular family. First, this is per- 
haps the largest family of graphs for which the Cauchy-Riemann operator admits 
a nice discretization. Indeed, restrictions of holomorphic functions to such graphs 
are preholomorphic to higher orders. This was the reason for the introduction of 
complex analysis on rhombic lattices by Duffin Duf68 in late sixties. More re- 
cently, the complex analysis on such graphs was studied for the sake of probabilistic 
applications |Mer01l IKen02l ICS08] . 

On the other hand, this seems to be the largest family where certain lattice 
models, including the Ising model, have nice integrability properties. In particular, 
the critical point can be defined with weights depending only on the local structure, 
and the star-triangle relation works out nicely. It seems that the first appearance 
of related family of graphs in the probabilistic context was in the work of Baxter 
|Bax78] . where the eight vertex and Ising models were considered on Z-invariant 
graphs, arising from planar line arrangements. These graphs are topologically 
the same as the isoradial ones, and though they are embedded differently into 
the plane, by [KS05] they always admit isoradial embeddings. In |Bax78j Baxter 
was not passing to the scaling limit, and so the actual choice of embedding was 
immaterial for his results. However, his choice of weights in the models would 
suggest an isoradial embedding, and the Ising model was so considered by Mercat 
|Mer01j . Boutilier and de Tiliere |BdT081 lBdT09] . Chelkak and the author [CS09] . 
Additionally, the dimer and the uniform spanning tree models on such graphs also 
have nice properties, see e.g. |Ken02) . 

We would also like to remark that rhombic lattices form a rather large family 
of graphs. While not every topological quadrangulation (graph all of whose faces 
are quadrangles) admits a rhombic embedding, Kenyon and Schlenker KS05 gave 
a simple topological condition necessary and sufficient for its existence. 

So this seems to be the most general family of graphs appropriate for our 
subject, and most of what we discuss below generalizes to it (though for simplicity 
we speak of the square and hexagonal lattices only). 



3. Applications of preholomorphic functions 

Besides being interesting in themselves, preholomorphic functions found several 
diverse applications in combinatorics, analysis, geometry, probability and physics. 
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After the original work of Kirchhoff, the first notable application was perhaps 
the famous article T3SST40 of Brooks, Smith, Stone and Tutte, who used preholo- 
morphic functions to construct tilings of rectangles by squares. 

Several applications to analysis followed, starting with a new proof of the Rie- 
mann uniformization theorem by Ferrand |LF55| . Solving the discrete version of 
the usual minimization problem, it is immediate to establish the existence of the 
minimizer and its properties, and then one shows that it has a scaling limit, which 
is the desired uniformization. Duffin and his co-authors found a number of similar 
applications, including construction of the Bergman kernel by Dieter and Mastin 
[DM71. There were also studies of discrete versions of the multi-dimensional com- 
plex analysis, see e.g. Kiselman's |KisQ5j . 

In |Thu86j Thurston proposed circle packings as another discretization of com- 
plex analysis. They found some beautiful applications, including yet another proof 
of the Riemann uniformization theorem by Rodin and Sullivan }RS87] . More in- 
terestingly, they were used by He and Schramm [HS93] in the best result so far 
on the Koebe uniformization conjecture, stating that any domain can be confor- 
mally uniformized to a domain bounded by circles and points. In particular, they 
established the conjecture for domains with countably many boundary compo- 
nents. More about circle packings can be learned form Stephenson's book |Ste05j . 
Note that unlike the discretizations discussed above, the circle packings lead to 
non-linear versions of the Cauchy- Riemann equations, see e.g. the discussion in 
|BMS05j . 

There are other interesting applications to geometry, analysis, combinatorics, 
probability, and we refer the interested reader to the expositions by Lovasz |Lov04j , 
Stephenson |Ste05j . Mercat }Mer07j . Bobenko and Suris jBS08j . 

In this note we are interested in applications to probability and statistical 
physics. Already the Kirchhoff's paper |Kir47] makes connection between the 
Uniform Spanning Tree and preharmonic (and so preholomorphic) functions. 

Connection of Random Walk to preharmonic functions was certainly known 
to many researchers in early twentieth century, and figured implicitly in many 
papers. It is explicitly discussed by Courant, Fricdrichs and Lewy in CFL28 , 
with preharmonic functions appearing as Green's functions and exit probabilities 
for the Random Walk. 

More recently, Kenyon found preholomorphic functions in the dimer model (and 
in the Uniform Spanning Tree in a way different from the original considerations 
of Kirchhoff). He was able to obtain many beautiful results about statistics of the 
dimer tilings, and in particular, showed that those have a conformally invariant 
scaling limit, described by the Gaussian Free Field, see |KenOO[ IKenOlj . More 
about Kenyon's results can be found in his expositions Ken04, Ken09]. An ap- 
proximately preholomorphic function was found by the author in the critical site 
percolation on the triangular lattice, allowing to prove the Cardy's formula for 
crossing probabilities [SmiOlb, Smi Olaj . 

Finally, we remark that various other discrete relations were observed in many 
integrable two dimensional models of statistical physics, but usually no explicit 
connection was made with complex analysis, and no scaling limit was considered. 
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Here we are interested in applications of integrability parallel to that for the Ran- 
dom Walk and the dimer model above. Namely, once a preholomorphic function is 
observed in some probabilistic model, we can pass to the scaling limit, obtaining 
a holomorphic function. Thus, the preholomorphic observable is approximately 
equal to the limiting holomorphic function, providing some knowledge about the 
model at hand. Below we discuss applications of this philosophy, starting with the 
Ising model. 



4. The Ising model 

In this Section we discuss some of the ways how preholomorphic functions appear 
in the Ising model at criticality. The observable below was proposed in [Smi06 for 
the hexagonal lattice, along with a possible generalization to O(N) model. Similar 
objects appeared earlier in Kadanoff and Ceva |KC71j and in Mercat |Mer01j . 
though boundary values and conformal covariance, which are central to us, were 
never discussed. 

The scaling limit and properties of our observable on isoradial graphs were 
worked out by Chelkak and the author in CSQ9]. It is more appropriate to consider 
it as a fermion or a spinor, by writing F(z)V~dz, and with more general setup one 
has to proceed in this way. 

Earlier we constructed a similar fermion for the random cluster representation 
of the Ising model, see Smi06, SmilO and our joint work with Chelkak CSQ9] for 
generalization to isoradial graphs (and also independent work of Riva and Cardy 
[RC06 for its physical connections). It has a simpler probabilistic interpretation 
than the fermion in the spin representation, as it can be written as the probability 
of the interface between two marked boundary points passing through a point 
inside, corrected by a complex weight depending on the winding. 

The fermion for the spin representation is more difficult to construct. Below we 
describe it in terms of contour collections with distinguished points. Alternatively 
it corresponds to the partition function of the Ising model with a y/z monodromy 
at a given edge, corrected by a complex weight; or to a product of order and 
disorder operators at neighboring site and dual site. 

We will consider the Ising model on the mesh e square lattice. Let f2 e be 
a discretization of some bounded domain O C C. The Ising model on Q e has 
configurations a which assign ±1 (or simply ±) spins o(v) to vertices v € Sl e and 
Hamiltonian defined (in the absence of an external magnetic field) by 

H(a) = - J2 <K«M«) . 

{u,v) 

where the sum is taken over all edges (it, v) inside f2 c . Then the partition function 
is given by 
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Figure 2. Left: configuration of spins in the Ising model with Dobrushin boundary 
conditions, its contour representation, and an interface between two boundary points. 
Right: an example of a configuration considered for the Fermionic observable: a number 
of loops and a contour connecting a to 2. It can be represented as a spin configuration 
with a monodromy at z. 



and probability of a given spin configuration becomes 

P(cr) = exp(-PH(a))/Z . 

Here /3 > is the temperature parameter (behaving like the reciprocal of the 
actual temperature), and Kramers and Wannier have established [KW41 that its 
critical value is given by j3 c = log + l) /2. 

Now represent the spin configurations graphically by a collection of interfaces - 
contours on the dual lattice, separating plus spins from minus spins, the so-called 
low-temperature expansion, see Figure [5] A contour collection is a set of edges, 
such that an even number emanates from every vertex. In such case the contours 
can be represented as a union of loops (possibly in a non-unique way, but we do not 
distinguish between different representations). Note that each contour collection 
corresponds to two spin collections which are negatives of each other, or to one if 
we fix the spin value at some vertex. The partition function of the Ising model can 
be rewritten in terms of the contour configurations u as 

g ^ y ^length of contours 

Each neighboring pair of opposite spins contributes an edge to the contours, and 
so a factor of x = exp(— 2/3) to the partition function. Note that the critical value 
is x c = exp(— 2/3 c ) = y/2 — 1. 

We now want to define a preholomorphic observable. To this effect we need 
to distinguish at least one point (so that the domain has a non-trivial conformal 
modulus). One of the possible applications lies in relating interfaces to Schramm's 
SLE curves, in the simplest setup running between two boundary points. To obtain 
a discrete interface between two boundary points a and b, we introduce Dobrushin 
boundary conditions: + on one boundary arc and — on another, see Figure [5] 
Then those become unique points with an odd number of contour edges emanating 
from them. 
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windings W = 1 windings tt W = — i 




windings ±2tt W = — 1 winding— 37r W = i 




Figure 3. Examples of Fermionic weights one obtains depending on the winding of the 
interface. Note that in the bottom left example there are two ways to trace the interface 
from a to z without self-intersections, which give different windings ±2n, but the same 
complex weight W = — 1. 

Now to define our fermion, we allow the second endpoint of the interface to 
move inside the domain. Namely, take an edge center z inside f2 £ , and define 

F e (z) := J2 xlCnSth ° f contom ' s ^ z)) , (7) 

w(a— ^z) 

where the sum is taken over all contour configurations ui = uj(a — > z) which have 
two exceptional points: a on the boundary and z inside. So the contour collection 
can be represented (perhaps non-uniquely) as a collection of loops plus an interface 
between a and z. 

Furthermore, the sum is corrected by a Fermionic complex weight, depending 
on the configuration: 

W(uj(a — > z)) := exp (— i s winding (7, a — > z)) . 

Here the winding is the total turn of the interface 7 connecting a to z, counted in 
radians, and the spin s is equal to 1/2 (it should not be confused with the Ising 
spins ±1). For some collections the interface can be chosen in more than one way, 
and then we trace it by taking a left turn whenever an ambiguity arises. Another 
choice might lead to a different value of winding, but if the loops and the interface 
have no "transversal" self-intersections, then the difference will be a multiple of 4w 
and so the complex weight VV is well-defined. Equivalently we can write 

W(u(a -> z)) = X* signcd turns of 7 , A := exp , 

see Figure [3] for weights corresponding to different windings. 
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Figure 4. Fisher graph for a region of the square lattice, a spin configuration and a 
corresponding dimer configuration, with dimers represented by the bold edges. 

Remark 1 . Removing complex weight W one retrieves the correlation of spins on 
the dual lattice at the dual temperature x* , a corollary of the Kramers- Wannier 
duality. 

Remark 2. While such contour collections cannot be directly represented by 
spin configurations, one can obtain them by creating a disorder operator, i.e. a 
monodromy at z: when one goes one time around z, spins change their signs. 

Our first theorem is the following, which is proved for general isoradial graphs 
in |CS09j , with a shorter proof for the square lattice given in [CSlOj : 

Theorem 1 (Chelkak, Smirnov). For Ising model at criticality, F is a preholo- 
morphic solution of a Riemann boundary value problem. When mesh e — > 0, 



where P is the complex Poisson kernel at a: a conformal map Q — > C + such that 
a t— > oo. Here both sides should be normalized in the same chart around b. 

Remark 3. For non-critical values of x observable F becomes massive preholomor- 
phic, satisfying the discrete analogue of the massive Cauchy-Riemann equations: 
dF = im{x - x c )F, cf. |MS09j . 

Remark 4. Ising model can be represented as a dimer model on the Fisher graph. 
For example, on the square lattice, one first represents the spin configuration as 
above — by the collection of contours on the dual lattice, separating + and — 
spins. Then the dual lattice is modified with every vertex replaced by a "city" of 
six vertices, see Figure |U It is easy to see that there is a natural bijection between 



F e (z)/V~e =t VP^z) 



inside 
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contour configurations on the dual square lattice and dimer configuration on its 
Fisher graph. 

Then, similarly to the work of Kenyon for the square lattice, the coupling func- 
tion for the Fisher lattice will satisfy difference equations, which upon examination 
turn out to be another discretization of Cauchy-Riemann equations, with different 
projections of the preholomorphic function assigned to six vertices in a "city" . One 
can then reinterpret the coupling function in terms of the Ising model, and this is 
the approach taken by Boutilier and de Tiliere |BdT08| lBdT09] . 

This is also how the author found the observable discussed in this Section, 
observing jointly with Kenyon in 2002 that it has the potential to imply the con- 
vergence of the interfaces to the Schramm's SLE curve. 

The key to establishing Theorem Q] is the observation that the function F 
is preholomorphic. Moreover, it turns out that F satisfies a stronger form of 
preholomorphicity, which implies the usual one, but is better adapted to fermions. 

Consider the function F on the centers of edges. We say that F is strongly 
(or spin) preholomorphic if for every centers u and v of two neighboring edges 
emanating from a vertex w, we have 

Proj(F(u), 1/V5) = Proj(F(u), 1/y/a) , 

where a is the unit bisector of the angle uwv, and Proj(p, q) denotes the orthogonal 
projection of the vector p on the vector q. Equivalently we can write 

F(v) +aF(v) = F{u) + aF{u) . (8) 

This definition implies the classical one for the square lattice, and it also easily 
adapts to the isoradial graphs. Note that for convenience we assume that the 
interface starts from a in the positive real direction as in Figure [21 which slightly 
changes weights compared to the convention in [CS09 . 

The strong preholomorphicity of the Ising model fermion is proved by con- 
structing a bijection between configurations included into F(v) and F(u). Indeed, 
erasing or adding half-edges wu and wv gives a bijection uj o uj between config- 
uration collections {uj(u)} and {uj(v)}, as illustrated in Figure [5j To check ©, 
it is sufficient to check that the sum of contributions from uj and uj satisfies it. 
Several possible configurations can be found, but essentially all boil down to the 
two illustrated in Figure \5\ 

Plugging the contributions from Figure [5] into the equation ((5J , we are left to 
check the following two identities: 

A + AA = 1 + Al , Xx + X~Xx = A 2 + AA 2 . (9) 

The first identity always holds, while the second one is easy to verify when x = 
x c = y2 — 1 and A = cxp(— iri/A). Note that in our setup on the square lattice 
A (or the spin s) is already fixed by the requirement that the complex weight is 
well-defined, and so the second equation in (j9]) uniquely fixes the allowed value of 
x. In the next Section we will discuss a more general setup, allowing for different 
values of the spin, corresponding to other lattice models. 
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contributes A C\ to F(v) 



k: 



contributes 1 C\ to F(u) 



contributes Xx C2 to F(v) 




contributes A 2 C2 to F(u) 



v I y a 
I w u ) 



Figure 5. Involution on the Ising model configurations, which adds or erases half-edges 
vw and uw. There are more pairs, but their relative contributions are always easy to 
calculate and each pair taken together satisfies the discrete Cauchy-Riemann equations. 
Note that with the chosen orientation constants Ci and C2 above are real. 



To determine F using its preholomorphicity, we need to understand its behavior 
on the boundary. When z G <9f2 e , the winding of the interface connecting a to z 
inside il £ is uniquely determined, and coincides with the winding of the boundary 
itself. This amounts to knowing Arg(_F) on the boundary, which would be sufficient 
to determine F knowing the singularity at a or the normalization at b. 

In the continuous setting the condition obtained is equivalent to the Riemann 
Boundary Value Problem (a homogeneous version of the Riemann- Hilbert-Privalov 
BVP) 



Im (F(z) ■ (tangent to dQ) 



1/2 







(10) 



with the square root appearing because of the Fermionic weight. Note that the 
homogeneous BVP above has conformally covariant solutions (as %/dz-forms) , and 
so is well defined even in domains with fractal boundaries. The Riemann BVP 
(TTU|) is clearly solved by the function \J P^(z), where P is the Schwarz kernel at a 
(the complex version of the Poisson kernel), i.e. a conformal map 



p ■. n 



a <-> 00 



Showing that on the lattice F e satisfies a discretization of the Riemann BVP 
(jlOP and converges to its continuous counterpart is highly non-trivial and a pri- 
ori not guaranteed - there exist "logical" discretizations of the Boundary Value 
Problems, whose solutions have degenerate or no scaling limits. We establish con- 
vergence in [CS09 by considering the primitive f* F 2 (u)du, which satisfies the 
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Dirichlct BVP even in the discrete setting. The big technical problem is that in 
the discrete case F 2 is no longer preholomorphic, so its primitive is a priori not 
preholomorphic or even well-defined. Fortunately, in our setting the imaginary 
part is still well-defined, so we can set 



1 r z 

H e (z) := — Im / F{zfdz . 



While the function H is not exactly preharmonic, it is approximately so, vanishes 
exactly on the boundary, and is positive inside the domain. This allows to complete 
the (at times quite involved) proof. A number of non-trivial discrete estimates is 
called for, and the situation is especially difficult for general isoradial graphs. We 
provide the needed tools in a separate paper [CS08 . 

Though Theorem Q] establishes convergence of but one observable, the latter 
(when normalized at b) is well behaved with respect to the interface traced from 
a. So it can be used to establish the following, see jCSlOj : 

Corollary 1. As mesh of the lattice tends to zero, the critical Ising interface in 
the discretization of the domain with Dobrushin boundary conditions converges 
to the Schramm's SLE(3) curve. 

Convergence is almost immediate in the topology of (probability measures on 
the space of) Loewner driving functions, but upgrading to convergence of curves 
requires extra estimates, cf. [KS09, DCHN09] ICS10] . Once interfaces are related 
to SLE curves, many more properties can be established, including values of di- 
mensions and scaling exponents. 

But even without appealing to SLE, one can use preholomorphic functions to a 
stronger effect. In a joint paper with Hongler [HS10] we study a similar observable, 
when both ends of the interface are allowed to be inside the domain. It turns out 
to be preholomorphic in both variables, except for the diagonal, and so its scaling 
limit can be identified with the Green's function solving the Riemann BVP. On 
the other hand, when two arguments are taken to be nearby, one retrieves the 
probability of an edge being present in the contour representation, or that the 
nearby spins are different. This allows to establish conformal invariance of the 
energy field in the scaling limit: 

Theorem 2 (Hongler, Smirnov). Let a £ Q and (x € ,y e ) be the closest edge from 
a G fi e . Then, as e — > 0, we have 

a 2 7T 

w reel V2 l n (a) 

Efree [° x °y\ = ^ — ' £ + O (e) , 

where the subscripts + and free denote the boundary conditions and Iq is the ele- 
ment of the hyperbolic metric on fl. 

This confirms the Conformal Field Theory predictions and, as far as we know, 
for the first time provides the multiplicative constant in front of the hyperbolic 
metric. 
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Figure 6. The high-temperature expansion of the O(N) model leads to a gas of disjoint 
simple loops. Probability of a configuration is proportional to AT* loopB i length . We study 
it with Dobrushin boundary conditions: besides loops, there is an interface between two 
boundary points a and 6. 

These techniques were taken further by Hongler in [HonlOj . where he showed 
that the (discrete) energy field in the critical Ising model on the square lattice 
has a conformally covariant scaling limit, which can be then identified with the 
corresponding Conformal Field Theory. This was accomplished by showing con- 
vergence of the discrete energy correlations in domains with a variety of boundary 
conditions to their continuous counterparts; the resulting limits are conformally 
covariant and are determined exactly. Similar result was obtained for the scaling 
limit of the spin field on the domain boundary. 



5. The O(N) model 

The Ising preholomorphic function was introduced in |Smi06] in the setting of 
general O(N) models on the hexagonal lattice. It can be further generalized to a 
variety of lattice models, see the work of Cardy, Ikhlef, Rajabpour |RC07[ ITC09j . 
Unfortunately, the observable seems only partially preholomorphic (satisfying only 
some of the Cauchy-Riemann equations) except for the Ising case. One can make 
an analogy with divergence-free vector fields, which are not a priori curl-free. 

The argument in the previous Section was adapted to the Ising case, and some 
properties remain hidden behind the notion of the strong holomorphicity. Below 
we present its version generalized to the O(N) model, following our joint work 
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[DCSlOj with Duminil-Copin. While for N ^ 1 we only prove that our observable 
is divergence- free, it still turns out to be enough to deduce some global information, 
establishing the Nienhuis conjecture on the exact value of the connective constant 
for the hexagonal lattice: 

Theorem 3 (Duminil-Copin, Smirnov). On the hexagonal lattice the number C(k) 
of distinct simple length k curves from the origin satisfies 

lim ylogC(fc) = log J 2 + V2 . (11) 

k— >oo k 

Self-avoiding walks on a lattice (those without self-intersections) were proposed 
by chemist Flory [Flo53 as a model for polymer chains, and turned out to be an 
interesting and extensively studied object, see the monograph MS93 . 

Using Coulomb gas formalism, physicist Nienhuis argued that the connective 
constant of the hexagonal lattice is equal to v2 + \/2, meaning that (ITTj) holds. 
He even proposed better description of the asymptotic behavior: 

C(k) « f^ + y/^j /c 11/32 , k^oo. (12) 

Note that while the exponential term with the connectivity constant is lattice- 
dependent, the power law correction is supposed to be universal. 

Our proof is partially motivated by Nienhuis' arguments, and also starts with 
considering the self-avoiding walk as a special case of O(N) model at N = 0. While 
a "half-preholomorphic" observable we construct does not seem sufficient to imply 
conformal invariance in the scaling limit, it can be used to establish the critical 
temperature, which gives the connective constant. 

The general O(N) model is defined for positive integer values of N, and is a 
generalization of the Ising model (to which it specializes for N = 1), with ±1 spins 
replaced by points on a sphere in the ./V-dimensional space. We work with the 
graphical representation, which is obtained using the high-temperature expansion, 
and makes the model well defined for all non-negative values of N. 

We concentrate on the hexagonal lattice in part because it is trivalent and so 
at most one contour can pass through a vertex, creating no ambiguities. This 
simplifies the reasoning, though general graphs can also be addressed by introduc- 
ing additional weights for multiple visits of vertices. We consider configurations 
uj of disjoint simple loops on the mesh e hexagonal lattice inside domain Q e , and 
two parameters: loop-weight N > and (temperature-like) edge-weight x > 0. 
Partition function is then given by 

g ^ y jy# *°°P S j^ cn &ii of contours 

A typical configuration is pictured in Figure |6l where we introduced Dobrushin 
boundary conditions: besides loops, there is an interface 7 joining two fixed bound- 
ary points a and b. It was conjectured by Kager and Nienhuis [KN04 that in 
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the interval N £ [0,2] the model has conformally invariant scaling limits for 
x = x c (N) := 1/V2 + V2 - N and x £ (x c (N),+oo). The two different limits 
correspond to dilute/dense regimes, with the interface 7 conjecturally converging 
to the Schramm's SLE curves for an appropriate value of n £ [8/3, 4] and K £ [4, 8] 
correspondingly. The scaling limit for low temperatures x £ (0, x c ) is not confor- 
mally invariant. 

Note that for N — 1 we do not count the loops, thus obtaining the low- 
temperature expansion of the Ising model on the dual triangular lattice. In partic- 
ular, the critical Ising corresponds to x = 1/ V3 by the work |Wan50] of Wannier, 
in agreement with Nienhuis predictions. And for x — 1 one obtains the critical site 
percolation on triangular lattice (or equivalently the Ising model at infinite temper- 
ature). The latter is conformally invariant in the scaling limit by |Smi01b|ISmi01a] . 

Note also that the Dobrushin boundary conditions make the model well-defined 
for N = 0: then we have only one interface, and no loops. In the dilute regime 
this model is expected to be in the universality class of the self-avoiding walk. 

Analogously to the Ising case, we define an observable (which is now a para- 
fermion of fractional spin) by moving one of the ends of the interface inside the 
domain. Namely, for an edge center z we set 

F e (z) := xlCUgth ° f contours W(uj(a -> z)) , (13) 

uj{a— ¥z) 

where the sum is taken over all configurations uj — ui{a — > z) which have disjoint 
simple contours: a number of loops and an interface 7 joining two exceptional 
points, a on the boundary and z inside. As before, the sum is corrected by a 
complex weight with the spin seK: 

W(u)(a —> 2)) := exp (— i s winding(7, a — > z)) , 
equivalently we can write 

W(w(a->«)) = A # signcd turns of 7 , A := exp ■ 

Note that on hexagonal lattice one turn corresponds to 7r/3, hence the difference 
in the definition of A. 

Our key observation is the following 

Lemma 4. For N £ [0, 2], set 2 cos (9) — N with parameter 6 £ [0, tt/2] . Then for 

x- l = 2 cos = ^2 - \/2^1v : , or (14) 

x- 1 = 2cos(^) = V2 + \/2^1v : , (15) 

the observable F satisfies the following relation for every vertex v inside fi e : 

(p - v)F(p) + (q - v)F(q) + (r - v)F(r) = 0, (16) 
where p, q, r are the mid-edges of the three edges adjacent to v. 



4ir 

7T + 36> 
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Figure 7. To obtain the parafermionic observable in the O(N) model we consider config- 
urations with an interface joining a boundary point z to an interior point z and weight 
them by a complex weight depending on the winding of the interface. 

Above solution (fT4|) corresponds to the dense, and (|15|) - to the dilute regime. 
Note that identity (fTo| is a form of the first Kirchhoff 's law, but apart from the 
Ising case N = 1 we cannot verify the second one. 

To prove LemmaHl we note that configurations with an interface arriving at p, q 
or r can be grouped in triplets, so that three configurations differ only in immediate 
vicinity of v, see Figure [5] It is enough then to check that contributions of three 
configurations to (|16[) sum up to zero. But the relative weights of configurations 
in a triplet are easy to write down as shown in Figure [51 and the coefficients in the 
identity (1161) are proportional to the three cube roots of unity: 1, r := cxp(i27r/3), 
f (if the neighbors of v are taken in the counterclockwise order). Therefore we 
have to check just two identities: 

iV + TA 4 +fA 4 = 0, 
l + r2;A + fxA = 0. 

Recalling that A = exp (—isw/3), the equations above can be recast as 

± (tt - 6) + 2irk , k G Z , 



2tt 



4s- 
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Ci to F(p) xX Ci to F{q) xX d to F{r) 

a 



N C 2 to F{p) A 4 C 2 to F(q) A 4 C 2 to F(r) 




Figure 8. Configurations with the interface ending at one of the three neighbors of v are 
grouped into triplets by adding or removing half-edges around v. Two essential examples 
of triplets are pictured above, along with their relative contributions to the identity (|13[) . 

The first equation implies that 

± (-I + 5)-5-I*- (17 » 

and the second equation then determines the allowed value of x uniquely. Most of 
the solutions of (|17|) lead to observables symmetric to the two main ones, which 
are provided by solutions to the equations (fT4|) and f|15|> . 

When we set N = 0, there are no loops, and configurations contain just an 
interface from a to z, weighted by x lcngth . This corresponds to taking 6 = ir/2 

and one of the solutions is given by s = 5/8 and x c = l/vHvl, as predicted 
by Nienhuis. To prove his prediction, we observe that summing the identity (|16[) 
over all interior vertices implies that 

where the sum taken over the centers z of oriented edges rj(z) emanating from the 
discrete domain f2 e into its exterior. Since F(a) = 1 by definition, we conclude 
that F for other boundary points sums up to 1. As in the Ising model, the winding 
on the boundary is uniquely determined, and (for this particular critical value of 
x), one observes that considering the real part of F we can get rid of the complex 
weights, replacing them by explicit positive constants (depending on the slope of 
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the boundary). Thus we obtain an equation 

^ ^ ^ ^ ^length of contours 

z£dCl e \{a} u(a— tz) 



regardless of the size of the domain f2 e . A simple counting argument then shows 
that the series 



E c w xk = E 



length 

o 

simple walks from a inside C 



converges when x < x c and diverges when x > x c , clearly implying the conjecture. 

Note that establishing the holomorphicity of our observable in the scaling limit 
would allow to relate self-avoiding walk to the Schramm's SLE with k = 8/3 and 
together with the work LSW04 of Lawler, Schramm and Werner to establish the 
more precise form (|12|) of the Nienhuis prediction. 



6. What's next 

Below we present a list of open questions. As before, we do not aim for complete- 
ness, rather we highlight a few directions we find particularly intriguing. 

Question 1. As was discussed, discrete complex analysis is well developed for 
isoradial graphs (or rhombic lattices) , see }Duf68| IMerOlj IKen02[ ICS08] . Is there 
a more general discrete setup where one can get similar estimates, in particular 
convergence of preholomorphic functions to the holomorphic ones in the scaling 
limit? Since not every topological quadrangulation admits a rhombic embedding 
[KS05], can we always find another embedding with a sufficiently nice version of 
discrete complex analysis? Same question can be posed for triangulations, with 
variations of the first definition by Isaacs ([3]) , like the ones in the work of Dynnikov 
and Novikov [DN03] being promising candidates. 

Question 2. Variants of the Ising observable were used by Hongler and Kytola 
to connect interfaces in domains with more general boundary conditions to more 
advanced variants of SLE curves, see [HK09] . Can one use some version of this ob- 
servable to describe the spin Ising loop soup by a collection of branching interfaces, 
which converge to a branching SLE tree in the scaling limit? Similar argument 
os possible for the random cluster representation of the Ising model, see :KS10 . 
Can one construct the energy field more explicitly than in [HonlOj . e.g. in the 
distributional sense? Can one construct other Ising fields? 

Question 3. So far "half-preholomorphic" parafermions similar to ones discussed 
in this paper have been found in a number of models, see jSmi06, RC06, RC07, 
IIC09j . but they seem fully preholomorphic only in the Ising case. Can we find the 
other half of the Cauchy-Riemann equations, perhaps for some modified definition? 
Note that it seems unlikely that one can establish conformal invariance of the 
scaling limit operating with only half of the Cauchy-Riemann equations, since 
there is no conformal structure present. 
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Question 4. In the case of the self-avoiding walk, an observable satisfying only 
a half of the Cauchy-Riemann equations turned out to be enough to derive the 
value of the connectivity constant IDCSlOj . Since similar observables are available 
for all other O(N) models, can we use them to establish the critical temperature 
values predicted by Nienhuis? Our proof cannot be directly transfered, since some 
counting estimates use the absence of loops. Similar question can be asked for 
other models. 

Question 5. If we cannot establish the preholomorphicity of our observables 
exactly, can we try to establish it approximately? With appropriate estimates 
that would allow to obtain holomorphic functions in the scaling limit and hence 
prove conformal invariance of the models concerned. Note that such more gen- 
eral approach worked for the critical site percolation on the triangular lattice 
[SmiOlbl ISmiOla) . though approximate preholomorphicity was a consequence of 
exact identities for quantities similar to discrete derivatives. 

Question 6. Can we find other preholomorphic observables besides ones men- 
tioned here and in [Smi06 ? It is also peculiar that all the models where preholo- 
morphic observables were found so far (the dimer model, the uniform spanning 
tree, the Ising model, percolation, etc.) can be represented as dimer models. Are 
there any models in other universality classes, admitting a dimer representation? 
Can then Kenyon's techniques |Ken04[ IKen09j be used to find preholomorphic 
observables by considering the Kasteleyn's matrix and the coupling function? 

Question 7. Throughout this paper we were concerned with linear discretizations 
of the Cauchy-Riemann equations. Those seem more natural in the probabilistic 
context, in particular they might be easier to relate to the SLE martingales, cf. 
|Smi06] . However there are also well-known non-linear versions of the Cauchy- 
Riemann equations. For example, the following version of the Hirota equation for 
a complex-valued function F arises in the context of the circle packings, see e.g. 
|BMS05j : 

(F(z + ie) - Fjz - e)) (F(z - ie) - F(z + e)) = 

(F(z + ie)-F(z + e))(F(z-ie)-F(z-e)) ' [ ' 

Can we observe this or a similar equation in the probabilistic context and use 
it to establish conformal invariance of some model? Note that plugging into the 
equation (|18|) a smooth function, we conclude that to satisfy it approximately it 
must obey the identity 

{d x F(z)f + (d v F{z)f = . 

So in the scaling limit (|18|) can be factored into the Cauchy-Riemann equations 
and their complex conjugate, thus being in some sense linear. It does not seem 
possible to obtain "essential" non-linearity using just four points, but using five 
points one can create one, as in the next question. 

Question 8. A number of non-linear identities was discovered for the correlation 
functions in the Ising model, starting with the work of Groeneveld, Boel and 
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Kasteleyn GBK78] IBK78] . We do not want to analyze the extensive literature 
to-date, but rather pose a question: can any of these relations be used to define 
discrete complex structures and pass to the scaling limit? In two of the early 
papers by McCoy, Wu and Perk }MW80| IPer80j , a quadratic difference relation 
was observed in the full plane Ising model first on the square lattice, and then 
on a general graph. To better adapt to our setup, we rephrase this relation for 
the correlation C(z) of two spins (one at the origin and another at z) in the Ising 
model at criticality on the mesh e square lattice. In the full plane, one has 

C{z + ie)C(z - ie) + C{z + e)C(z - e) = 2C(z) 2 . (19) 

Note that C is a real-valued function, and the equation (1191) is a discrete form of 
the identity 

C(z)AC(z) + \VC(z)\ 2 = . 

The latter is conformally invariant, and is solved by moduli of analytic functions. 
Can one write an analogous to (|19p identity in domains with boundary, perhaps 
approximately? Can one deduce conformally invariant scaling limit of the spin 
correlations in that way? 

Question 9. Recently there was a surge of interest in random planar graphs and 
their scaling limits, see e.g. [DS09, LGPQ8]. Can one find observables on random 
planar graphs (weighted by the partition function of some lattice model) which 
after an appropriate embedding (e.g. via a circle packing or a piecewise-linear 
Riemann surface) are preholomorphic? This would help to show that planar maps 
converge to the Liouville Quantum Gravity in the scaling limit. 

Question 10. Approach to the two-dimensional integrable models described here 
is in several aspects similar to the older approaches based on the Yang-Baxter 
relations |Bax89 . Some similarities are discussed in Cardy's paper |Car09j . Can 
one find a direct link between the two approaches? It would also be interesting to 
find a link to the three-dimensional consistency relations as discussed in BMS09 . 

Question 11. Recently Kenyon investigated the Laplacian on the vector bundles 
over graphs in relation to the spanning trees [KenlO . Similar setup seems natural 
for the Ising observable we discuss. Can one obtain more information about the 
Ising and other models by studying difference operators on vector bundles over the 
corresponding graphs? 

Question 12. Can anything similar be done for the three-dimensional models? 
While preholomorphic functions do not exist here, preharmonic vector fields are 
well-defined and appear naturally for the Uniform Spanning Tree and the Loop 
Erased Random Walk. To what extent can they be used? Can one find any other 
difference equations in three-dimensional lattice models? 
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